%sample code
%example linking Java with matlab
%function PatchTestT3
clc
clear
close all
%have to call init_mesh2d first
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% create mesh
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
l = 1;
h = 1;
vertex = [0 0; l 0;l h;0 h];
edge = [1 2;2 3;3 4;4 1];
hdata.hmax = 0.25;
[node,element] = mesh2d(vertex,edge,hdata);
ms = createMesh(node,element);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% create element
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
M = fem2.Model(ms);
E = 2;
nu = 0.3;
ss = fem2.StressState.PLANE_STRESS;
cm = fem2.StVenantKirchhoffCM2D(E, nu, ss);
for i=1:ms.countFaces();
    e = fem2.T3Element(ms.getFace(i-1),cm);
    M.addElement(e);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% set boundary condition
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
ncx = find(node(:,1)==0);
ncy = find(node(:,2)==0);
nco = find(node(:,1)==0 & node(:,2)==0);
c1 = fem2.Constraint([0 0]);
c2 = fem2.Constraint([0 1]);
c3 = fem2.Constraint([1 0]);
for i=1:size(ncx(:),1);
    ms.getNode(ncx(i)-1).setConstraint(c2);
end
for i=1:size(ncy(:),1);
    ms.getNode(ncy(i)-1).setConstraint(c3);
end
for i=1:size(nco(:),1);
    ms.getNode(nco(i)-1).setConstraint(c1);
end

fcx = find(node(:,1)==l);
F = 0.5;
for i=1:size(fcx(:),1)-1;
    n1 = ms.getNode(fcx(i)-1);
    n2 = ms.getNode(fcx(i+1)-1);
    e = ms.addEdge(fem2.Edge([n1 n2]));
    ldx = math2.LinearFunctionRToR(-1,1,F,F);
    ldy = math2.LinearFunctionRToR(-1,1,0,0);
    el = fem2.EdgeLoad(e,[ldx ldy]);
    M.addLoad(el);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% solve
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
M.preAnalysis();
n = ms.countDofs();
kk = inf.jlinalg.Array2DMatrix(n,n);
M.assembleKt(kk);
r = zeros(1,n);
r = M.assembleR(r);
u = kk.getValues()\r;
M.setX(u);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% print results
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
node_deform = node;
for i=1:ms.countNodes();
    du = ms.getNode(i-1).getUHatLocal();
    node_deform(i,:) = node_deform(i,:)+du';
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% plot results
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% figure
% plot(node(:,1),node(:,2),'b*');
% hold on
% for i=1:size(node,1)
%     text(node(i,1),node(i,2),num2str(i));
% end
% plot(node_deform(:,1),node_deform(:,2),'r*');
plotDeformedMesh(M.getMesh());
